Power-law spin correlations in a perturbed honeycomb spin model 
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We consider spin-^ model on the honeycomb lattice [l[ in presence of a weak magnetic field 
h a <C 1. Such a perturbation destroys exact integrability of the model in terms of gapless fermions 
and static Z2 fluxes. We show that it results in appearance of a long-range tail in the irreducible 
dynamic spin correlation function: ((s z (t, r)s z (0, 0))) oc h 2 f(t,r), where f(t,r) tx [max(t, r)]~ 4 is 
proportional to the density polarization function of fermions. 



Quantum spin liquids, QSL's (see e.g. Refs. 0-0]) 
present examples of strongly correlated quantum phases 
which do not develop any kind of local order, while their 
specific entropy vanishes at zero temperature. Critical, 
or algebraic QSL's are characterized by spin correlation 
functions that decay as some power of distance and time. 
In some cases, the correlation asymptotics can be de- 
duced from a representation of spin operators in terms of 
almost-free fermions [f|. However, a complete calculation 
based on a microscopic Hamiltonian has not been demon- 
strated due to the lack of suitable exactly solvable models 
(in more than one spatial dimension). We show in the 
present Letter that the anisotropic spin- ^ model on the 
honeycomb lattice, proposed by one of us [lj , can be used 
as a starting point for the construction of an analytically 
treatable critical QSL. This result may seem surprising 
since it is known Q that the original model [l| possesses 
no spin correlations at the distances longer than a single 
lattice bond. We will see however, that a small perturba- 
tion of the model [H, e.g. weak external magnetic field, 
is sufficient to "turn on" long-range spin correlations, al- 
beit with a small overall prefactor. Thus we disagree with 
the statement made in Ref. Q that short-range charac- 
ter of spin correlations survives in the presence of a weak 
magnetic field. 

We consider the model defined by the Hamiltonian: 

rl = J ^2 (trjO-i) (a-jtii) - y^hjO-j. (1) 

l=(ij) i 

Unit vectors n; are parallel to x, y and z axis for the 
corresponding links x, y and z of the honeycomb lattice, 
see Fig.1. At = the Hamiltonian ([1]) was solved ex- 
actly |l| via a mapping to a free fermion Hamiltonian. 
In this approach, each spin cr,; is represented in terms 
of four Majorana operators a, cf , cf , cf with the fol- 
lowing anticommutation relations: {cf,c^} = 28^5^, 
so that <rf = icicf. In terms of these new operators, 
the zero-field Hamiltonian reads rl — — c i u ij c j 
and Uij = icfc" are constants of motion: [H, ity] = 0, 
with = ±1. The ground state \G) corresponds to 
an arbitrary choice of {u^} that minimizes the energy. 
It is convenient to introduce the notion of Z2 flux, de- 
fined for each hexagon 7r as a product <f>^ = Y[ Uij (since 




FIG. 1: A honeycomb lattice fragment with the z-z link that 
belongs to a given elementary cell. The indicated vector qjr = 
x shows the direction of oscillations found in the (s z (0)s z (r)} 
correlation function. 



Uij = —Uji, we have to choose a particular ordering in 
this definition: i £ even sublattice, j € odd sublattice). 
The ground state of this model is a symmetrized sum 
of states with different sets of integrals of motion {u^}, 
corresponding to all fluxes equal to 1. For practical cal- 
culations of physical quantities, one does not have to im- 
plement such symmetrization and can assume that all 
Uij = 1. We denote by H the corresponding Majorana 
Hamiltonian: H = —iJYluj) c i c j- It can be diagonalized 
with the use of Fourier transformation. The spectrum of 
the resulting free fermions is gapless and has two conic 
points. To begin, we recall the calculation of the spin- 
spin correlation function = ^crf (t)cr^(0)\ in the un- 
perturbed model with h = Q • Spin operator af acting 
on the ground state produces two Z2 fluxes and creates 
a fermion. Since states with different flux configurations 
are mutually orthogonal and fluxes do not move in the 
process of time evolution governed by the Hamiltonian 
rl, a non-zero result for the correlation function is possi- 
ble only if the second spin operator er^ creates the same 
pair of fluxes. Thus the sites i and j are either the same 
or nearest neighbors. For larger separations between i 
and j, one has gij = 0. However, this result is due to 
the static nature of Z2 fluxes; furthermore, the ground 
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state is a linear combination of states with the same flux 
pattern. Both these features are destroyed by any pertur- 
bation (for example, external magnetic field) that does 
not commute with operators {iiij}. 

In this paper, we consider the honeycomb lattice model 
with an external magnetic field, which is treated as a 
weak perturbation. Before delving into calculations, we 
note that a generic magnetic field opens a gap A in the 
fermionic spectrum with A ~ h x h y h z / J 2 . In what 
follows we neglect this gap. This is definitely possible if 
one of the field components vanishes, i.e. if the field is 
directed in one of the coordinate planes. For the generic 
field direction, our results for spin correlations are ap- 
plicable for intermediate distances, 1 « r « 1/A. For 
simplicity we discuss the total spin in the r-th elementary 
cell, s r = a r i + o>,2 an d calculate z—z correlations only. 
Since the external magnetic field induces finite magne- 
tization, (s z ) 7^ 0, we study the irreducible correlation 
function: g(t,r) = ((sj;(t)sg(0))). We are interested in 
the long-time and/or long-distance asymptotics of g(t, r). 

It is convenient to introduce complex bond fcrmions, 
defined as follows: ip r = | (c r ,i + ic r , 2 ) and 4> r = 
\ (c z 1 + ic z 2 ) . Operator 4> r creates two fluxes in the 
plaquettes adjacent to the z-link in the elementary cell 
r. Note, that </>^0 r = 1+ ^' J , and hence the ground-state 
wave function \G) satisfies (j)^4 l r \G) = \G). Spin opera- 
tor factorizes in the following way: s z = 2itp"(j>", where 
vectors tp" = (?/Vj ipt ) ; 4% = {4'ri't'r) are introduced 
and summation over a is implied. In the absence of mag- 
netic field, the correlation function of flux operators <p r is 
local: G<£ (r, t) = (0+(£)(/> o (O)) = ip(t)8 r fi, which leads to 
locality of the spin correlations. However, once the mag- 
netic field is turned on, one finds G^ (r, t) to be nonzero 
and proportional to h 2 at any r, which leads to spin cor- 
relation at large distances. We start from the expression 
for (s^(i)sg(O)), expanded up to the second order in h z : 

(s z r (t)4m = -y E JJ dTldT2 ' (2) 
■{Ts; (t)4 (0)< (nK 2 (r 3 )>. 

The irreducible correlation function equals g (t, r) = 
(*r(*)«o(0)> - <*o (0)> 2 , where 

(s z o (0)) = -ih z J^ I dT(T s ~(0) s *(T)). (3) 

r •* 

Thus we have to calculate two-spin and four-spin cor- 
relation functions. For these correlation functions to be 
non-zero, the flux configuration which results from the 
action of the two (four) spin operators on the ground 
state should coincide with the original flux configuration. 
In particular, the two-spin correlator vanishes unless 
r = 0, so that we have (Tsq (0) s z r (r)) = u (r) 5 r $ with 
u(t) = (Ts^ (0) Sq (t)). It is clear that u(— t) = u(t). 



The expression for the magnetization thus simplifies: 

/>oo 

(s z ) = -2ih z / u(t)cLt. (4) 
Jo 

Similarly, 

(Ts z (t)s z (0)< (ri)< (r a )) = 
= h (r,t,Ti,T 2 ) 5 rur 5 r2 .o + /2 (r,t, tx, t 2 ) 5 ru oS r2 , r , 

where 

fi(r, t,rx,r 2 ) = (Ts z (t) s z (0) s z (n) s z (r 2 )) , (5) 
h(r, t, tx,t 2 ) - (Ts z (t) s z (0) s z (tx) s z (r a )) . 

In the t — > oo limit, the leading contributions to /i and 
f 2 come from the regions t\ m t, t 2 and t\ ~ 0, 
t-2 t, respectively. We will see that the product of spin 
operators at nearby times, e.g. s z (t)s z (T2) in the second 
case, reduces to the product of two fermion operators 
(up to some renormalization). It follows that the four- 
spin correlation function is asymptotically proportional 
to the density polarization function of free fcrmions. 

Note that the Wick theorem is not directly applicable 
to spin averages because each spin operator creates both 
a fermion and some flux, the latter acting as a scatter- 
ing potential for propagating fermions. To proceed with 
the calculation, one has to rewrite spin operators s z in 
equations ([5]) in terms of fermions ip and and then 
move the cj> operators to the right, commuting them with 
exponential evolution factors. To this end, we use the 
identities 

where the Hamiltonian H r differs from the original 
Hamiltonian H by inverting sign of the u variable which 
belongs to the z-link in the elementary cell r: H r = 
H + V r , where V r = 4 J (V'r Vv — !)• I n this way, all 
spin correlators can be represented as correlators of non- 
interacting fermions in the presence of external time- 
dependent potential. The calculation of u(t) is a sim- 
ple task discussed in Refs. 0, 0|- Using the identity 
e iHt e -iH r t _ y CX p V^(r)dr^, one arrives at the 

following result: 

«(t) = 4 (TiP (r) iP+ (0) e- l ti M T ') dT ' J> for r > 0. 

The next step is to calculate fx,2(r, t,Tx,T 2 ). We con- 
sider explicitly all different time ordcrings in the expres- 
sion ([5]); it is enough to choose t > 0, since g r (—t) = 

g*- r (ty- 

1) r 2 > tx > t > 0; 2) r 2 > t > > n; 

3) t > > t 2 > tx; 4) t > t 2 > tx > 0; 

5) t>T 2 >0>Tx; 6) t 2 > t > tx > 0, 
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while other 6 domains 1'..6' can be obtained by the per- for the particular time domain j = 2. We get for the 
mutation T\ < — > r 2 . following expression (summation over a..S is implied): 

Let us illustrate how to perform the calculation of \ 



n 



(2) 



= (*g (r 2 ) s z r (t) s g (0) si (n)) = 16 ( e ^Vo^oe- iffr2 e iH Vf^e- iff V^^ iHTl ^^e 
= 16 ^e i ^ T2 V'oe"^ oT2 e i - fft) *V're" < - ffr - Vje iH '- T1 V'J.e-^ T1 44^ ^ 
= -16 (r^ (r 2 ) W (t) V>o + (0) 1&+ (ti) e"^ 2 '^) . 



(7) 



To proceed from the first to the second line, we used 
the fact that the only relevant sequence of superscripts 
is a/3j5 = 1122 (recall that <p+<p r \G) = \G) , while 
<p r <p+ \G) = and & = <pt 2 = 0). Similarly, for / 2 (2) 
one obtains: 

/ 2 (2) = 16 (t^, (r a ) W i'o (0) (n) e"*/ ^W*) . 

(8) 

The potentials V^r) are piece-wisc-constant functions 
of time which can be easily read off the order of fcrmionic 
operators in ([7]), ©: 
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In the same way exact expressions for f? , analogous 
to (|7I8|) . can be obtained for all other time domains 1...6. 
However, they can hardly be evaluated exactly in the 
closed form. The problem of their calculation resembles 
the one encountered while exploring the Fermi Edge Sin- 
gularity problem f]~oj ] , so we can analyze it similarly. 

The representation of spin correlation functions in the 
form (J7J allows us to use the Wick theorem for fermions, 
which makes a diagrammatic expansion of fl over the 
potential V-~ possible. Note that apart from the nor- 
mal Green function G (t, r) = (Tip (r, t) ip + (0, 0)), the 
anomalous Green function F (t, r) = (Tip (r, t) ip (0, 0)) = 
(Tip + (0,0) ip + (r,t)) has also to be taken into account 
(we calculate both of them below). The sum of all dia- 
grams for each of is of the form = l&e C i • 
where the first factor is the the sum of closed-loop dia- 
grams, and the second factor uf is the a sum of open- 
line diagrams. 

The closed-loop contribution equals e c '> <J> ( t ' ri ' T2 ' = 
^Te~ l I V ' (3) ^ dT ^ . In the limit of large time separation 

between pairs of points {t, n}, {0, r 2 } or {t, r 2 }, {0, n}, 

(i) 

the asymptotic form of C\ can be simply determined. 



For example, the expressions for C\ 2 read: 

C[ 2) ps -i(Q(t + T 2 -T 1 ) + Sn r t), 
C| 2) ps -% CI (r 2 - n - t) . 

In this equation, fl is the energy of the fcrmionic ground 
state in the presence of two adjacent fluxes (Ct ps 0.04J, 
see [l[), while SQ r stands for the interaction energy of two 
flux pairs separated by distance r. Therefore, the factor 
cxp (cf^ (i, T\, r 2 )^j rapidly oscillates with frequency fi. 

Each term in the sum of open-line diagrams corre- 
sponds to a particular pairing of four fcrmionic operators 

_____ / 2 \ 

in the product ([7]) or ([5]). For example, L\ ' is given by 
the following equation: 

L { 2 2) = (TiP r (r 2 ) ip+ {t))f (Ti, Q (0) ip+ (n)) 2 2) - (9) 

- (TW (r 2 ) iP (0)) 2 2) (TiPt (i) iPt (rr)) 2 2) - 

- (TiP r (r a ) (n)) 2 2) (T^o (0) W + W)f ■ 

where (T...)^ stands for (r...e- i S y P* r \ . 

Finally, the spin correlation function g(t, r) is given by 
the time integral over n, t 2 of the oscillating function 
. The contribution from the domains 2 + 2' reads: 

/0 /* oo 
/ dT X dT 2 e c ^ hf . (10) 

Up to this point, all calculations have been exact for any 
r, t. To proceed further, we have to make some approx- 
imations. We use the inequality t ^> J~ x and average 

over fast oscillations of gf\t,v) as a function of t. It 

fa) 

is easy to see that g\ has a slowly varying part be- 

(2) 

cause (r) = for r G (0,i). On the other hand, 

(2) 

g[ is purely oscillating and vanishes upon the averag- 

f 2 ) 

ing. Considering the expression (|10p for <? 2 , we find 
that the main contribution to it comes from t\ ~ 
and r 2 ps t (the result of integration is determined by 
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a small neighborhood of the boundary points due to os- 
cillations of the integrand), so the corresponding expres- 
sion in ((HJ is of the form of fermionic density-density 

(2) 

correlation function. Now we have to calculate L\ ' for 
Ti«0 and r 2 ~ t. Note that for such time arguments, 

(2) 

the external potential V 2 as a function of r turns on 
for two short intervals (of the order of O -1 ), while the 
separation between the pulses is large, t ^> J -1 . In 
this case the long-time (or large-distance) asymptotics 

of the correlation function reads (Ttjj r (t 2 ) i/jq { t i)) 2 ~ 
G r (t) <f) (0 — Ti) <j) i T 2 — t) , where </>(r), 4>(t) are some di- 
mensionless functions of Jr. The double integral in (fTU)) 
is thus factorized, and the renormalization due to the 
functions (f>(t), <fi(t) adds an overall numerical coefficient 
only, which we denote by h$ 2 (it is the same for the con- 
tributions from all time domains). Calculating the dom- 
inant first term in ^ (the other oscillate as functions of 

(2) 

t), we obtain g 2 ■ 



these asymptotics to Eq. ((TTj) . we obtain the final result: 



^ (2) = -16/1= x 



f poo \ 2 

92 ' = -^ n z X [J u(t)<1t\ - 
-lQh 2 Jh 2 [F (t, r) F (-t, -r) + G (t, r) G (-t, -r)] 



Similar considerations are applicable for the other time 
domains 1..6 show that all relevant contributions have 
a similar feature: the integration over {ti,T2} is domi- 
nated by some neighborhood of points and t. Collecting 
everything and subtracting (sg (0)) , we obtain: 
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-f [F (t, r) F (-t, -r) + G (t, r) G (-t, -r)] 

flri 



(11) 

where free fermion Green functions G and F are cal- 
culated below in Eq. (fT2|) . Thus we have found that 
the spin correlation function g (t, r) is proportional to 
the density-density correlation function of band fcrmions, 
with the coefficient oc h 2 z . The parameter Hq in Eq. (TlTj) 
can be estimated (up to a numerical constant) to be 
h ~ J ■ 

Now we turn to the calculation of fermionic Green 
functions G(t,r) and F (t, r) which enter (fTTj) . The ex- 
pression for "vector" composed of these Green functions 
in the energy-coordinate representation reads: 



(G (r) F (r)) = — ^ (( £ + 5 /p ) cos(pr), -5R/ P sin(pr)) 



I/pI 



iS 



(12) 

where /(p) = 2iJ (l + e' pn i + e lpn 2) and n lj2 = 
^±i,^^ in the standard (x,y) coordinates. We ex- 
pand /(p) near the conical point K = ^§7r, -^71"^ to get 
long-time behavior of G (t, r) , and F (t, r) . Substituting 



0(f,r) 



1G 

Z2 



hz \ 2 (r 2 -3(Jtf) cos 2 (fq^r) 



(13) 

where = x is the unit vector along x, orthogonal to 
the direction of z — z link. The singularity of this expres- 
sion at r/t = \/3 J is cut off by the finite width of the 
Brillouin zone, which was sent to infinity while calculat- 
ing integrals leading to Eqs. (fT3"| . The above result may 
seem surprising due to the apparent anisotropy demon- 
strated by fast oscillations in expression (fl~3|) as a func- 
tion of the x-component of r. This anisotropy is due to 
our choice to calculate correlations of the z components 
of the spin in the unit cell. Similar calculation for x or 
y components lead to analogous results with anisotropy 
vectors q^r y , which are perpendicular to the correspond- 
ing lattice links. The correlation function (fl~3|) was calcu- 
lated in the lowest nontrivial order over perturbation h z . 
To account for higher-order terms in this expansion, one 
should study the effects of flux motion upon the fermion 
polarization function. 

In conclusions, we have shown that under weak per- 
turbation due to magnetic field, spin operators acquire 
nonzero projection ~ h 2 on the density of band fcrmions, 
thus long-range spin correlations appear. Therefore 
weakly perturbed honeycomb spin model may be consid- 
ered as an example of the critical QSL. We expect that 
the same mechanism of the coupling of spins to fermion 
density can be realized for similar models on the deco- 
rated honeycomb lattice [1, Q . For the gapful state [§] , 
we expect spin correlations to decay with a correlation 
length/time determined by the gap in the fermion spec- 
trum, whereas for the spin metal state 0] asymptotic 
behavior f(t,r) cx [max(t,r)p is expected. Another 
perturbation leading to the similar effect is a weak modi- 
fication of the triad of n vectors, which will also produce 
nonzero spin correlations at arbitrary distances. 
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